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Abstract 

Background: Protein phosphorylation is one of the most important post-translational modifications involved in 
the regulation of plant growth and development as well as diverse stress response. As a member of the Poaceae, 
Brachypodium distachyon L. is a new model plant for wheat and barley as well as several potential biofuel grasses 
such as switchgrass. Vegetative growth is vital for biomass accumulation of plants, but knowledge regarding the role of 
protein phosphorylation modification during vegetative growth, especially in biofuel plants, is far from comprehensive. 

Results: In this study, we carried out the first large-scale phosphoproteome analysis of seedling leaves in Brachypodium 
accession Bd21 using Ti0 2 microcolumns combined with liquid chromatography-tandem mass spectrometry 
(LC-MS/MS) and MaxQuant software. A total of 1470 phosphorylation sites in 950 phosphoproteins were identified, 
and these phosphoproteins were implicated in various molecular functions and basic cellular processes by gene ontology 
(GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analyses. Among the 950 phosphoproteins 
identified, 127 contained 3 to 8 phosphorylation sites. Conservation analysis showed that 93.4% of the 950 
phosphoproteins had phosphorylation orthologs in other plant species. Motif-X analysis of the phosphorylation 
sites identified 13 significantly enriched phosphorylation motifs, of which 3 were novel phosphorylation motifs. 
Meanwhile, there were 91 phosphoproteins with both multiple phosphorylation sites and multiple phosphorylation 
motifs. In addition, we identified 58 phosphorylated transcription factors across 21 families and found out 6 significantly 
over-represented transcription factor families (C3H, Trihelix, CAMTA, TALE, MYB_related and CPP). Eighty-four protein 
kinases (PKs), 8 protein phosphatases (PPs) and 6 CESAs were recognized as phosphoproteins. 

Conclusions: Through a large-scale bioinformatics analysis of the phosphorylation data in seedling leaves, a 
complicated PKs- and PPs- centered network related to rapid vegetative growth was deciphered in B. distachyon. 
We revealed a MAPK cascade network that might play the crucial roles during the phosphorylation signal transduction 
in leaf growth and development. The phosphoproteins and phosphosites identified from our study expanded our 
knowledge of protein phosphorylation modification in plants, especially in monocots. 
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Background 

Protein phosphorylation is one of the common reversible 
post-translational modifications involved in the regulation 
of plant growth and diverse processes [1], Various tech- 
niques have been developed for specific enrichment of 
phosphopeptides, among which Fe + -IMAC (Immobilized 
metal affinity chromatography) [2,3], Ti0 2 -MOAC (Metal 
oxide affinity chromatography) [4,5] and SCX (Strong 
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cation exchange) chromatography [6,7] are the three most 
common used methods for large-scale phosphoproteome 
studies. In recent years, through these phosphopeptide 
enrichment techniques combined with high-accuracy 
mass spectrometry (MS) and related bioinformatics, 
many large-scale phosphoproteomic analyses were per- 
formed in different plant species, such as Ambidopsis 
thaliana [8-11], Oryza sativa [10], Medicago truncatula 
[12,13], Glycine max [11], Zea mays [14,15] and Triticum 
aestivum [16,17]. However, knowledge of protein phos- 
phorylation modification in other plant species, especially 
for members of the Poaceae, is far from comprehensive. 
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Vegetative growth rapidly increases the photosyn- 
thetic capacity and size of plants, and is vital for 
biomass accumulation, especially in biofuel plants. 
Brachypodium distachyon L., a member of the Pooi- 
deae subfamily and a temperate wild annual grass 
endemic to the Mediterranean and Middle East [18], 
has rapidly become a model plant, especially for poten- 
tial biofuel grasses such as switchgrass {Panicum virga- 
tum L.). It possesses many attractive attributes such as 
a small diploid genome of 272 Mbp, short growth 
cycle, self-fertility and simple nutrient requirements 
[19] as well as competence to be efficiently trans- 
formed [20]. 

Leaf expansion is a major aspect of plant vegetative 
growth. It increases light capture, which powers photo- 
synthesis and thus biomass production. Considerable 
works have investigated the mechanisms of leaf expan- 
sion [21-24], but little is known about post-translational 
phosphorylation modification of proteins during leaf 
expansion in the period of rapid vegetative growth. 
During this period the leaf cell wall undergoes dynamic 
changes to allow cells to expand, but at the same time 
cells must maintain the mechanical strength required to 
resist the forces of turgor pressure [25]. During plant 
growth and development, leaf cell numbers and size are 
rapidly increasing and cell walls must adapt to these 
changes. The components and mechanisms of under- 
lying signaling systems to achieve this process remain 
largely unknown, but emerging evidences have impli- 
cated several receptor-like kinases as regulators of cell 
wall function [26-30]. 

Interactions of protein kinases (PKs) with their sub- 
strates are to a large extent determined by residues sur- 
rounding the phosphorylation sites and the pattern of 
residues is named as a phosphorylation motif [31,32]. 
Motif-X analysis of large phosphorylation site data sets 
can detect significantly enriched phosphorylation motifs 
and predict the corresponding PKs [33]. Schwartz et al. 
[34] used Motif-X to determine phosphorylation motifs 
in yeast, fly, mouse, and man. Trost et al. [35] found 
eight specific phagosomal phosphorylation motifs in- 
duced by IFN-y. Bennetzen et al. [36] used Motif-X to 
analyze phosphorylation sites during DNA damage re- 
sponses and identified a novel [sxxQ] motif, confirming 
that Motif-X is a powerful program for revealing signifi- 
cant and novel phosphorylation motifs. 

Recently, a large-scale proteomic and phosphoproteo- 
mic study of B. distachyon leaves under salt stress was 
performed in our laboratory [37]. In the present study, 
we focused on phosphorylation modification during 
seedling leaf growth and development at the omics level 
for the first time in B. distachyon and revealed a com- 
plicated phosphorylation signal transduction network 
during rapid leaf growth and development. 



Results and discussion 

Identification of phosphorylation sites and phosphoproteins 

In this study, a large-scale phosphoproteome analysis in 
B. distachyon seedling leaves was performed to explore 
the complex protein phosphorylation network of signal- 
ing and regulatory events. The strategy used in this 
study is shown in Figure 1. To increase the number of 
identified phosphopeptides and eliminate false positives, 
three biological replicates were used for phosphopeptide 
enrichment and LC-MS/MS analysis. A total of 1937 
phosphopeptides containing 2449 phosphorylation sites 
were identified (Additional file 1: Table SI). The raw 
mass spectrometry proteomics data have been deposited 
to the ProteomeXchange Consortium (http://proteome- 
central.proteomexchange.org) via the PRIDE partner 
repository [38] with the dataset identifier PXD000868. 
Only the phosphopeptides with unambiguous (class I) 
phosphorylation sites identified from three biological 
repeats were used for further analysis. Finally, the 1470 
unambiguous phosphorylation sites distributed in 1367 
phosphopeptides, corresponding to 950 phosphoproteins, 
were screened and used for further analysis (Additional 
file 1: Table SI, sheet A). Among the 1470 unambiguous 
phosphorylation sites, 1313 Ser, 155 Thr, and 2 Tyr phos- 
phorylation sites were identified, accounting for 89.32%, 
10.54%, and 0.14% of the total sites, respectively (Figure 2, 
panel A). For the lower identification of phosphorylated 
tyrosine, there are three possible reasons: first, we didn't 
perform focused analysis during MS runs in this study, 
which may lead to the reduced discovery of the low abun- 
dant Tyr phosphopeptides; second, the recovery rate of 
phospho-Tyr peptides could be different because of the 
enrichment techniques and data processing softwares; 
third, different organs, development stages or stress treat- 
ments can also lead to the variation of the phospho-Tyr 
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Figure 1 Strategy for a large-scale phosphoproteomics study 
on seedling leaves of B. distachyon. 
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Figure 2 Analysis of identified phosphorylated sites. (A) Distribution of the phosphorylated amino acids. (B) Distribution of phosphorylation 
sites in the phosphopeptides. (C) Distribution of phosphorylation sites in the phosphoproteins. 



identification. The review by Ghelis [39] also gave similar 
explanations for the lower phospho-Tyr proportion in 
plants. The distributions of phosphorylation sites in phos- 
phopeptides and phosphoproteins were shown in Figure 2, 
panel B and C, respectively. Among the 1367 phospho- 
peptides, 1266 (92.61%) of them each contained one 
phosphorylation sites, only 99 (7.24%) contained two 
and 2 (0.15%) contained three phosphorylation sites 
(Figure 2, panel B). From the phosphoprotein perspective, 
656 (69.05%) and 167 (17.58%) of the 950 phosphopro- 
teins were detected to have 1 and 2 phosphorylation sites, 
respectively. One hundred and twenty seven phospho- 
proteins contained more than two phosphorylation sites, 
of which 25 had at least five phosphorylation sites. Two 
phosphoproteins (Bradilg66870.1 and Bradi3g53907.1) 
each possessed 8 phosphorylation sites (Figure 2, panel C). 
The representative mass spectra of all the seven phospho- 
peptides of Bradilg66870.1 containing the 8 phosphory- 
lated sites were shown in Additional file 2: Figure SI. 

In order to obtain an overview of phosphorylation 
events in seedling leaves of B. distachyon, all of the 950 
identified phosphoproteins were used to perform GO, 
Pfam annotation. Finally, 840 of them were annotated by 
Blast2GO and 804 of them possessed Pfam domains 
(Additional file 1: Table SI, sheet B). Significantly (FDR 
adjusted p < 0.05) enriched GO items of biological process, 
molecular function and cellular component are shown in 
Figure 3. Both the percent and p-value of each item were 
displayed to evaluate these GO items. From the biological 



process perspective, "nucleobase, nucleoside, nucleotide 
and nucleic acid metabolic process (23.05%, GO:0006139, 
FDR:2.6E - 43)", "response to stimulus (19.37%, GO:0 
050896, FDR:1.5E - 106)", "Protein modification process 
(15.36%, GO:0006464, FDR:3.8E - 17)", "cellular compo- 
nent organization (14.84%, GO:0016043, FDR:2.3E - 98)", 
"transport (14.21%, GO:0006810, FDR:4.5E - 23)", and 
"cellular amino acid and derivative metabolic process 
(11.26%, GO:0006519, FDR9.9E - 66)" were over-represen 
ted (Figure 3, panel A). From the cellular component per- 
spective, "membrane (32.63%, GO:0016020, FDR:8.7E - 
99)" and "nucleus (9.37%, GO:0005634, FDR:1.5E - 8)" 
were significantly enriched (Figure 3, panel B). From 
the molecular function perspective, "protein binding 
(25.57%, GO:0005515, FDR:6.1E - 86)", "nucleotide bind- 
ing (24.42%, GO:0000166, FDR3.1E - 22)", "kinase activity 
(11.58%, GO:0016301, FDR:1.3E - 8)" and "DNA binding 
(11.37%, GO:0003677, FDR:3E-12)" were the highly 
enriched molecular function items (Figure 3, panel C). 

Conservation and KEGG pathway analyses of 
phosphoproteins 

To investigate the conservation of phosphoproteins 
between B. distachyon and other plant species, the 
sequences of the 950 B. distachyon phosphoproteins 
were used as queries to blast a phosphoprotein data- 
base which were constructed based on datasets from 
P3DB (Plant Protein Phosphorylation DataBase) [40], 
MORE (Medicago-Omics Repository) [13] and PhosPhAt 
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Figure 3 Gene ontology enrichment analysis of the 950 identified phosphoproteins. (A) Biological processes. (B) Cellular components. 
(C) Molecular functions. 



4.0 [41]. Oryza sativa and Arabidopsis thaliana, as the 
model species of monocots and dicots, respectively, pos- 
sess more comprehensive phosphorylaton information 
than other species. Recently, a large phosphoprotein data- 
set of the model legume Medicago truncatula was released 
[13]. B. distachyon was therefore compared with the three 
species to determine the degree of conservation of the 
phosphoproteins among them. The threshold was set as 
Score > 80, E-value < IE- 10 and Identity > 30%. Finally, 
582 (61.3%) phosphoproteins were identified with phos- 
phorylation orthologs in all the three species, 202 (21.3%) 
had phosphorylation orthologs in two of the three species, 
and 103 (10.8%) had phosphorylation orthologs in only 
one species (Additional file 3: Table S2). Only 63 phos- 
phoproteins had no phosphorylation orthologs in the 
three species. The conserved analysis results are shown 
in Additional file 3: Table S2. Further biological process 
enrichment analysis of the 582 highly conserved phospho- 
proteins (Additional file 4: Figure S2, panel A) showed 
that "cellular protein modification process (GO:000 
6464, FDR:2.4E-2)" and "signal transduction (GO:0007165, 
FDR:5.0E-3)" were significantly over-represented from 
the total identified phosphoproteins (Additional file 4: 
Figure S2, panel B). 

To reveal the phosphoprotein-associated pathways, 
the 950 phosphoproteins were mapped to KEGG; 31 
B. distachyon KEGG pathways with not less than 3 hits 



were highlighted by the phosphoproteins identified in 
this study (Additional file 5: Table S3 and Additional 
file 6: Figure S3). Many fundamental biological pathways 
were highlighted by the mapping of phosphoproteins 
from this work, including carbohydrate metabolism in- 
volving carbon metabolism (Additional file 6: Figure S3, 
panel D), glycolysis/gluconeogenesis (panel G), pyruvate 
metabolism (panel N) and pentose phosphate pathway 
(panel P), energy metabolism including carbon fixation 
(panel I), starch and sucrose metabolism (panel L) and 
oxidative phosphorylation (panel R), and nucleotide me- 
tabolism including purine metabolism (panel M). For 
genetic information processing, transcription in spliceo- 
somes (panel A), translation related to RNA transport 
(panel C), ribosome (panel E) and mRNA surveillance 
(panel F) were highly enriched. Two signal transduction 
pathways were particularly evident: the plant hormone 
signal transduction pathway (panel Q) and the phos- 
phatidylinositol signaling system (panel AC). Further 
analysis of the 582 highly conserved phosphoproteins- 
related pathways showed that carbohydrate metabolism 
and energy metabolism-associated pathways, splicing, 
RNA transport and mRNA surveillance pathways stood 
out (Additional file 5: Table S3). The result of KEGG 
pathway analysis was consistent with that from GO ana- 
lysis, indicating that the methods used here were effect- 
ive and reliable. 
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Phosphorylation of transcription factors 

In growth and development of plants, many transcrip- 
tion factors (TFs) are affected through phosphorylation 
by PKs or dephosphorylation by protein phosphatases 
(PPs), which in turn either positively or negatively regulate 
TF activity in facilitating a program of gene expression 
that results in changed cell behavior [42]. Based on GO 
and KEGG pathway analyses, transcription-associated 
proteins were one of the most significant groups of all 
phosphoproteins, and TFs particularly were the main 
constituents. Fifty eight TFs belonging to 21 different 
families were identified (Table 1) based on alignments of 
our phosphoproteins with B. distachyon TFs in PlantTFDB 
(http://planttfdb.cbi.pku.edu.cn/index.php?sp=Bdi) [43]. In 
PlantTFDB, there are 1751 B. distachyon TFs classified into 
56 families. To determine the TF families significantiy regu- 
lated by phosphorylation, a TF family enrichment test was 
performed from 21 TF families contained 58 TFs compared 
with the background set (1751 TFs from 56 TF families 
in B. distachyon). The number of phosphorylated TFs 
belonging to a TF family was compared against the total 
number of TFs in that family of B. distachyon. Finally, 6 
significantly over-represented (p < 0.05) TF families 
(C3H, Trihelix, CAMTA, TALE, MYB_related, and CPP) 
were found (Figure 4, panel A). Bd21 is diploid inbred line 
with its genome size about 272 Mb and contains 26,552 
gene loci, coding for 31,029 distinct mRNA molecules 
[18]. Genes of all the 58 TFs were localized in the five 
Brachypodium chromosomes (Figure 4, panel B). The 
distribution of the six significantly over-represented TF 
families had no obvious chromosomal preference. 

Phosphorylation motif enrichment analysis 

To date, about 14,000 and 12,000 phosphorylation sites 
representing nearly 4100 and 4900 phosphoproteins in 
P3DB have been found in A. thaliana and O. sativa, 
respectively [40]. For A. thaliana, PhosPhAt 4.0 contains 
12,613 experimental phosphorylation sites and 5663 
phosphoproteins [41]. In addition, a large dataset for the 
model legume M. truncatula containing 13,506 phos- 
phorylation sites localized in 3926 phosphoproteins was 
available in MORE [13]. Phosphorylation modification of 
these sites are regulated by PKs, which account for about 
5.5% of the Arabidopsis genome [44]. Different PKs have 
preference for specific substrates, and many kinase- 
associated phosphorylation motifs have been deciphered 
[45-49]. The Motif-X online tool was used in this study 
to identify the phosphorylation motifs, and 13 different 
motifs were enriched, including 11 Ser and 2 Thr motifs 
(Figure 5 and Additional file 7: Table S4, sheet A). The 
phosphoproteins containing each of the 13 motifs (listed 
in Additional file 7: Table S4, sheet B) were further used in 
GO enrichment analysis. As shown in Additional file 8: 
Figure S4, DNA binding or RNA binding was enriched 



from phosphoproteins with the motifs 1, 2, 3, 6, 10 or 13, 
whereas translation factor activity was enriched from 
phosphoproteins with the motifs 9 or 10. Phosphoproteins 
containing motifs 1, 4 or 9 were mainly involved in the cell 
cycle and located in ribosomes. Lipid binding and receptor 
activity were over-represented from motifs 2 and 8, re- 
spectively, and tropism was enriched from phosphopro- 
teins containing the motif 12. Organic cyclic compound 
binding and heterocyclic compound binding functions 
were enriched from the novel motif 5. No GO terms were 
significantly enriched from phosphoproteins containing 
motifs 7 or 11. 

According to the literature and databases [45-49], 10 
of these motifs were specific to known PKs and 3 had no 
corresponding known kinases (Additional file 7: Table S4, 
sheet A). Ten known motifs included proline-directed mo- 
tifs such as [s/tP] (Motifs 6 and 13) and [Pxs/tP] (Motifs 1 
and 12), which were potential substrates of mitogen- 
activated PK (MAPK), cyclin-dependent kinase (CDK) and 
CDK-like kinase. Basic motifs [sPR] (Motif 2) and [sPxR] 
(Motif 4) were recognized by growth-associated histone 
kinase (GHK), CDK or cell division cycle 2 (CDC2) kinase. 
[Rxxs] (Motif 8) was recognized by calcium/calmodulin- 
dependent PK II (CaMK-II). [MxRxxs] (Motif 7) was rec- 
ognized by C-terminal Src Kinase-homologous Kinase 1 
(CHK1). Acidic motifs [sDxE] (Motif 9) and [sxD] (Motif 
10) were recognized by casein kinase-II (CK-II). The other 
3 phosphorylation motifs had no known specific kinases 
and were therefore regarded as novel phosphorylation 
motifs found in B. distachyon. Further analysis of these 
motifs showed that [RxxsP] (Motif 3), as a basic motif, 
is a combination of [sP] (Motif 6) and [Rxxs] (Motif 8); 
and this motif may be recognized by MAPK, CDK or 
CaMK. [GxxsP] (Motif 5) is similar to [GsP], which 
is recognized by glycogen synthase kinase 3 (GSK-3), 
MAPK or CDK5. [sxG] (Motif 11) may be recognized by 
CKII [50]. Further inspection of the 127 phosphoproteins 
containing multiple phosphorylation sites (Additional 
file 1: Table SI, sheet B) showed that 122 of them pos- 
sessed at least one of the 13 phosphorylation motifs 
(Additional file 7: Table S4, sheet C). Among them, 91 
phosphoproteins possessed two or more kinds of motifs. 
Phosphorylation sites of most phosphoproteins were 
motif- specific, but some other phosphoproteins each 
contained multiple phosphorylation sites belonging the 
same motif. Two phosphoproteins (Bradilg51560.1 and 
Bradi3g03690.1) each contained four phosphorylation 
sites specific to one motif (Additional file 7: Table S4, 
sheet C). This suggests that most phosphoproteins with 
multiple sites are regulated by multiple PKs, and a small 
part of them are activated by only one kind of PKs. The 
phosphoprotein Bradilg66870.1 contained eight phos- 
phorylation sites belonging to five kinds of motifs (Figure 6, 
panel A), which may be regulated by different PKs such as 
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Table 1 Phosphorylated transcription factors identified 
in this study 



Phosphoprotein 
no. 


Number of 

phosphorylation 

sites 


Phosphorylation 
sites 


TF family 


Bradi1g69480.1 


1 


S372 


ARR-B 


Bradi1g14510.1 


1 


S576 


B3 


Bradi2g57800.1 


1 


S161 


bHLH 


Bradi3g39927.1 


4 


S1 1 ;S30;S1 85; 
S258 


bHLH 


Bradi1g05480.1 


3 


S6;S71;S399 


bZIP 


Bradi3g07540.1 


1 


S32 


bZIP 


Bradi3g38200.1 


1 


S1 12 


bZIP 


Bradi3g41 980.1 


5 


S38;S114;S241; 
S257;S271 


bZIP 


Bradi4g41 890.1 


1 


S29 


bZIP 


Bradi3g51090.1 


1 


S73 


C2H2 


Bradi4g37540.1 


1 


S56 


C2H2 


Bradi5g1 8530.1 


1 


S81 


C2H2 


Bradilgl 1760.1 


4 


S208;S511;S524; 
S641 


C3H 


Bradi1g24270.1 


2 


S493;S563 


C3H 


Bradi1g32977.1 


1 


S353 


C3H 


Bradi2g09350.1 


1 


S435 


C3H 


Bradi2g38247.1 


2 


S397;S550 


C3H 


Bradi2g58770.1 


2 


S206;S286 


C3H 


Bradi3g04650.1 


1 


S258 


C3H 


Bradi3g09570.1 


1 


S38 


C3H 


Bradi4g06290.1 


2 


S437;S448 


C3H 


Bradi4g19010.1 


2 


S215;S219 


C3H 


Bradi5g02310.1 


1 


T35 


C3H 


Bradi5g10750.1 


5 


S19;S21;S38;S42;S46 


C3H 


Bradi1g21372.1 


3 


S329;S777;S929 


CAMTA 


Bradi1g71810.2 


2 


S594;S973 


CAMTA 


Bradi3g23800.1 


1 


S956 


CAMTA 


Bradi2g20347.1 


1 


S234 


CPP 


Bradi4g02410.1 


1 


S226 


CPP 


Bradi1g07630.1 


1 


S190 


G2-like 


Bradi3g20610.1 




S100;S358 


GeBP 


Bradi3g57800.1 


1 


T86 


GeBP 


Bradi2g46230.1 




S952 


HB-other 


Bradi4g35760.1 




S863 


HD-ZIP 


Bradi2g49860.1 




S230 


HSF 


Bradi4g07110.1 




S416 


MYB 


Bradi1g08470.1 




S204 


MYB_related 


Bradi1g351 90.1 




S265 


MYB_related 


Bradi2g1 8250.1 




S426 


MYB_related 


Bradi2g42870.1 




S105 


MYB_related 



Table 1 Phosphorylated transcription factors identified 



in this study (Continued) 


Bradi3g01 340.3 1 


S1081 


MYB_related 


Bradi5g07300.1 1 


S72 


MYB_related 


Bradi1g60030.1 3 


S4;S8;T1 24 


NF-YB 


Bradi4g38527.1 1 


S19 


NF-YB 


Bradi1g09615.1 1 


S421 


TALE 


Bradi1g12510.1 1 


S7 


TALE 


Bradi1g74060.1 1 


S448 


TALE 


Bradi4g25020.1 1 


S75 


TALE 


Dn^i^n £03 ID 1 1 

DraQogjyjzU. I 


J I JO 


I Lr 


Dld0lzy4OJZU. I Z 




Trihelix 


Bradi3g00697.1 3 


S282;S418;S665 


Trihelix 


Bradi3g46210.1 1 


S65 


Trihelix 


Bradi4g37730.1 1 


S224 


Trihelix 


BradiSgl 1070.1 1 


S78 


Trihelix 


Bradi1g23340.1 1 


S105 


WRKY 


Bradi3g 19640.1 4 


S332;T352;T354; 
S409 


WRKY 


Bradi3g52837.1 1 


S265 


ZF-HD 


Bradi4g30240.1 1 


S96 


ZF-HD 



MAPK, CDK and CK-II. Conservation analysis showed 
that phosphorylation orthologs in other species were also 
identified as phosphoproteins with many phosphosites, 
and 5 of the 8 identified phosphorylation sites in 
Bradilg66870.1 were also identified in other species whereas 
the other three were the newly identified phosphorylation 
sites present only in B. distachyon (Figure 6, panel B). 

Phosphorylation of PKs and PPs 

There are 989 PKs and 131 PPs in Arabidopsis according 
to PlantsP (http://plantsp.genomics.purdue.edu/) [51]. In 
this study, 84 PKs and 9 PPs were identified as phospho- 
proteins (Table 2). To identify the interactions between 
the phosphorylated PKs and PPs and their potential sub- 
strates, protein-protein interaction (PPI) analysis was con- 
ducted by using the Search Tool for Retrieval of Interacting 
Genes/Proteins (STRING). All 478 euKaryotic Orthologous 
Groups (KOGs) (Additional file 1: Table SI, sheet B) repre- 
senting 685 of the 950 phosphoproteins were used to con- 
struct the phosphoprotein interaction network. To improve 
the accuracy of PPI analysis, the confidence level (score) 
was set to a high value (0.900). Finally, a complicated PPI 
network centered by PKs and PPs was constructed and 
displayed by using Cytoscape software (Additional file 9: 
Figure S5). Depending on the homology comparison with 
the PKs and PPs of Arabidopsis in the PlantsP database, 
the phosphorylated PKs and PPs identified in our study 
were classified into six classes (Table 1): transmem- 
brane receptor kinases and related non-transmembrane 
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bHLH 
MYB 
bZIP 
C2H2 
WRKY 
MYB related 




Phosphorylated TFs 
Background 



25% 



B 



Scale 
(Mbp) 



20 - 



60 - 



Chr1 
74.8 



bZIP bradi1g05480 
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Figure 4 Analysis of phosphorylated transcription factors (TFs). (A) TF family enrichment analysis by using a two-sided test of the hypergeometric 
distribution. The number of phosphorylated TFs belonging to a TF family was compared against the total number of TFs in that family of B. distachyon. The 
six significantly over-represented TF families are marked with *p < 0.05 or **p < 0.01 . (B) Chromosomal localization of genes of all phosphorylated TFs on 
the five chromosomes of 6. distachyon. The chromosome numbers and their size (Mb) are indicated at the top of each bar. 



Lv et al. BMC Genomics 2014, 15:375 
http://www.biomedcentral.com/1471-2164/15/375 



Page 8 of 18 



1 

2 
3 
4 
5 
6 
7 




¥ ? * * t 



AAA L£ 

eep§y A 




8 ii 



¥ ? T * * 



T'f T 




y * t * r 



* V t * t 



mvm rlllil 




¥ ¥ t * * 




r y r i * * 



Figure 5 Significantly enriched phosphorylation motifs. Among these motifs, 10 are known phosphorylation motifs, and 3 are newly 
identified (enclosed into green boxes). Detailed information and putatively associated kinases are shown in Additional file 7: Table S4. 



kinases, ATNl/CTRl/EDRl/GmPK6-like kinases, casein 
kinase I, non-transmembrane PKs, other kinases, and 
plant phosphatases. 

MAPK cascades are universal signal transduction mod- 
ules in plants and are involved in responses to various 
biotic and abiotic stresses, hormones, cell division and 
developmental processes [52]. Activation of a MAPK cas- 
cade can lead to changes in gene expression relating to cell 
wall biogenesis [53,54]. MAPK cascades contained three 
PK families, the MAPK, MAPK kinase (MAP2K) and 
MAPK kinase kinase (MEKK) families. In this study, 
one MAPK (Bradilg49 100.1; KOG0660), one MAP2K 
(Bradilg51000.1; KOG0581) and five MEKK (Bradilg 
10970.1, Bradi4g22760.1, Bradi5gl0670.1, Bradi5gl81 
80.1, and Bradi5g24870.1; KOG0198) were identified as 
phosphoproteins, and conserved among B. distachyon, O. 
sativa, A. thaliana and M. truncatula (Additional file 3: 
Table S2), clearly indicating the fundamental roles of 
MAPK cascades in plants. A sub-network (Figure 7) 
centered by MEKK, MAP2K and MAPK was extracted 
from the whole PPI network mentioned above. From 
the sub-network, the MAPK family was the crucial PK 



in MAPK cascades and directly activated the substrates, 
so most phosphoproteins including many other kind of 
PKs (CDK, PHOT and STE20) interacted with MAPK. 
MAPK and many other PKs such as CK1 and SNT7 and 
PP2C were activated by MAP2K whereas MAP2K was 
activated by MEKK. Many phosphoproteins interacted 
simultaneously with MEKK and MAPK (Group 1) or 
MAP2K and MAPK (Group 2), and some interacted 
with MEKK, MAP2K and MAPK (Group 3). Several 
PKs and PPs, including GSK3, PDK1, DSPK and PPI, 
were also regulated by both MEKK and MAPK whereas 
only DSPP was regulated by both MAP2K and MAPK. 
Two phosphorylation sites (Table 2 and Additional 
file 1: Table SI, sheet A) were identified from BdMAPK 
(Bradilg49100.1) and one of them was phosphorylated 
at Tyr220. An Ambidopsis MAPK (AtMPK4) is phos- 
phorylated at a Tyr residue during the activation process 
[55]. A previous study revealed that plant MAPKs con- 
tained a conserved phosphorylation motif, [tD/Ey] [56]. 
Interestingly, the BdMAPK we identified also contains 
this motif. Besides, we identified phosphorylation of 
11 BdCDPKs and 4 BdCDKs, which may function in 
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Figure 6 Phosphorylation sites analysis of Bradi1g66870.1. (A) Domain structures of Bradi1g66870.1 and distribution of the 8 phosphorylated 
sites and related phosphorylation motifs on the sequence of Bradi1g66870.1. (B) Sequence alignments among Bradi1g66870.1 and its phosphorylation 
orthologs from 0. sativa, A. thaliana, M. truncotulo, and G. max. Phosphorylation sites identified in this study and based on the information from P3DB 
are highlighted in red. 



defense response and regulation of cell division and dif- 
ferentiation [57-59]. 

Phosphorylation of signal transduction proteins 

In plants, many signal pathways are involved in vegeta- 
tive growth. In the present study, KEGG pathway ana- 
lysis showed that the plant hormone signal transduction 
pathway (Additional file 6: Figure S3, panel Q) and 
phosphatidylinositol signaling system (Additional file 6: 
Figure S3, panel AC) played important roles through 
phosphorylation modification of some crucial proteins. 
The hormones include cytokinin, abscisic acid, ethylene 
and brassinosteroid which mainly regulate cell division, 
cell elongation and shoot initiation. A homolog of type-B 
Arabidopsis response regulator (B-ARR, Bradilg69480.1) 
was found in the cytokinin regulated pathway. B-ARR 
is a DNA-binding transcriptional regulator, whose activity 
is regulated by phosphorylation of its phospho-accepting 



receiver domain through the His-kinase-mediated His 
to Asp phosphorelay and lack of this domain disturbs 
the development of young seedlings [60]. Likewise, 
other conserved phosphoproteins including a bZIP TF 
(Bradi3g38200.1), two ethylene-insensitive protein (Bradi4g 
08380.1 and Bradi5g00700.1) and three serine/threonine 
PKs (Bradilg07620.1, Bradi3g32330.1 and Bradi2g56267.1) 
were involved in the abscisic acid, ethylene or brassi- 
nosteroid regulated signal transduction pathways. An- 
other signal transduction pathway shown by KEGG 
was the phosphatidylinositol signaling system (Additional 
file 6: Figure S3, panel AC) that is associated with cell 
proliferation and differentiation [61]. Phosphatidylino 
sitol-3,4,5-trisphosphate 3-phosphatase (PTEN, Bradi 
4g08080.1) acts as a lipid phosphatase to produce 
phosphatidylinositol 4,5-biphosphate (PtdIns(4,5)P2). Phos 
phatidylinositol-4- phosphate 5-kinase (Bradilg44047.1) 
is involved in PtdIns(4,5)P2 synthesis. Phosphoinositide 
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Table 2 Phosphorylated PKs and PPs identified in this study 



Phosphoprotein no. 



Phosphorylation site count 



Phosphorylation sites 



KOG no. 



Class 1 -transmembrane receptor kinases and related non-transmembrane kinases 

Bradi1g07010.1 1 

Bradi1g207S0.1 1 

Bradi1g35477.1 1 

Bradi1g58260.1 3 

Bradil g5921 0.1 1 

Bradi1g60397.2 2 

Bradil g72430.1 2 

Bradil g76660.1 1 

Bradi2g01180.1 1 

Bradi2g12120.1 1 

Bradi2g27350.1 1 

Bradi2g34630.1 1 

Bradi2g43l10.1 2 

Bradi2g47000.1 1 

Bradi3g08660.1 1 

Bradi3g32330.1 1 

Bradi3g40370.1 1 

Bradi3g44250.1 1 

Bradi3g49160.1 1 
Bradi3g53350.1 
Bradi4g35050.1 
Bradi5g03410.1 
Bradi5g1 2540.1 
Bradi5g13960.1 



2 

1 

3 
2 
1 

Class 2-ATN1/CTR1/EDR1/GmPK6 like kinases 

Bradil g00580.1 2 
Bradil g281 10.2 1 
Bradil g28950.2 3 
Bradil g30720.1 1 
Bradil g47570.1 1 
Bradi2g46340.1 2 
Bradi2g49790.1 1 
Bradi3g09170.1 3 
Bradi3g60210.1 1 
Bradi4g36880.1 2 
Bradi4g41 870.1 4 
Class 3-casein kinase I 
Bradi2g14310.1 1 
Bradi2g47850.1 1 
Bradi3g54860.1 1 
Bradi2g33037.2 1 
Class 4-non-transmembrane PKs 
Bradil g06300.1 1 



S372 
S735 
S697 

S741;S743;S733 
S837 

S345;S370 

S382;S389 

S770 

S133 

S149 

S444 

T263 

S654;T674 

S682 

S1031 

S233 

S93 

S508 

S693 

S546;S535 
S24 

S342;T666;T672 

S665;S32S 

S612 

S514;S503 
S12 

S1105;S761;S569 

S477 

T14 

S333;S410 
S129 

S478;S481;S400 
S811 

S27;S1 53 

S23;S37;S234;T417 

S415 
S417 
S456 
S10 

S137 



KOG 1187 
KOG1187 
KOG1187 
KOG1187 
KOG1187 
KOG1187 
KOG1187 
KOG1187 
KOG1187 
KOG1187 
KOG1187 
KOG1187 
KOG1187 
KOG1187 
KOG1187 
KOG1187 
KOG1187 
KOG1187 
KOG1187 
KOG1187 
KOG1187 
KOG1187 
KOG1187 
KOG1187 

KOG0192 
KOG0192 
KOG0192 
KOG0192 
KOG0192 
KOG0192 
KOG0192 
KOG0192 
KOG0192 
KOG0192 
KOG0192 

KOG 1163 
KOG 1163 
KOG 1163 
KOG 1164 

KOG0032 
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Table 2 Phosphorylated PKs and PPs identified in this study (Continued) 



Bradi1g24240.1 


1 


S524 


KOG0032 


Bradi1g35090.1 


4 


S102;S26;T36;T46 


KOG0032 


Bradi1g56970.1 


2 


S557;S85 


KOG0032 


Bradi2g1 5520.1 


, 


S537 


KOG0032 


Bradi4g24390.1 


1 


S504 


KOG0032 


Bradi4g39870.1 


] 


S538 


KOG0032 


Bradi4g40300.1 


! 


S509 


KOG0032 


Bradi5g1 9430.1 




S40;S549;S553;S16 


KOG0032 


Bradil gl 0970.1 


1 


S230 


KOG0198 


Bradi2g26547.1 


1 


S972 


KOG0198 


Bradi4g22760.1 


1 


S155 


KOG0198 


Bradi5g10670.1 


1 


S607 


KOG0198 


Bradi5g1 81 80.1 




S335;S337 


KOG0198 


Bradi5g24870.1 




S893 


KOG0198 


Bradil g26670.1 




S322 


KOG0201 


Bradil g42257.2 




S502;S336;S477 


KOG0581 


Bradil g51000.1 


1 


T25 


KOG0581 


Bradi1g75150.1 




T30 


KOG0581 


Bradil g77700.1 




S538;S562;S625 


KOG0582 


Bradi3g05890.1 




S330;T331 


KOG0582 


Bradi3g31 110.1 




S537 


KOG0582 


Bradil g07620.1 




S179;S183 


KOG0583 


Bradi2g56267.1 


1 


S158 


KOG0583 


Bradil g23970.1 


1 


S614 


KOG0584 


Bradi2g56460.1 


, 


S394 


KOG0592 


Bradil g1 8450.2 


1 


S375 


KOG0600 


Bradil g62440.1 




S77 


KOG0600 


Bradi2g1 2937.2 




S509 


KOG0600 


Bradi2g26680.1 


1 


T198 


KOG0600 


Bradi2g61 31 0.1 


2 


S93;T1 98 


KOG0600 


Bradil gl 1 380.1 


3 


S635;S654;T28 


KOG0606 


Bradi4g45310.1 


4 


S378;S385;S286;S297 


KOG0610 


Bradi5g07360.1 


3 


S508;S306;T304 


KOG0610 


Bradi3g10350.1 


2 


S466;S427 


KOG0614 


Bradi2g38590.1 


1 


Y235 


KOG0658 


Bradi2g26510.1 


2 


S337;T1 73 


KOG0659 


Bradil g491 00.1 


2 


T218;Y220 


KOG0660 


Bradi5g 1 3980.1 


3 


S125;S190;S693 


KOG0663 


Bradil g1 0530.1 


1 


S546 


KOG0667 


Bradi4g00230.1 


1 


T544 


KOG0670 


Bradi2g42320.1 


1 


S33 


KOG0671 


Class 5-other kinases 








Bradil g08540.1 


1 


S92 


K0G1151 


Bradil g26580.1 


1 


S30 


K0G1166 


Bradi2g1 7660.1 


3 


S545;T531;T540 


KOG0594 
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Table 2 Phosphorylated PKs and PPs identified in this study (Continued) 



Class 6-plant phosphatases 



Bradi2g36370.1 


2 


S847;S500 


KOG0374 


Bradi1g24400.1 


1 


S414 


KOG0698 


Bradi1g47710.1 


1 


S339 


KOG0698 


Bradi3g39540.1 


1 


S504 


KOG0698 


Bradi3g521 10.1 


2 


S121;S131 


KOG0700 


Bradi1g42810.1 


1 


S232 


KOG1716 


Bradi2g37450.1 


3 


S283;S46S;T92 


KOG1716 


Bradi4g08080.1 




S73;S113 


KOG2283 



phospholipase C (PLC, Bradilgl6560.1) mediates the 
production of second messenger molecules diacylglycerol 
(DAG) and inositol 1,4,5-trisphosphate (IP3) by hydrolysis 
of PtdIns(4,5)P2. The three enzymes were all identified as 
phosphoproteins (Additional file 1: Table SI). 

Phosphorylation of proteins related to cell wall expansion 

The Poaceae family can be distinguished from most other 
land plants by the composition of their cell walls [62,63]. 
B. distachyon is a good model for cell wall (mostly primary 
wall) studies of temperate grasses [24,64]. Cell wall loosen- 
ing is essential for leaf expansion during vegetative growth 
of plants, and this must be balanced with polymer synthe- 
sis and wall restrengthening to prevent the cell wall from 
rupture. Multiple receptor-like PKs (RLKs) are implicated 
in cell wall signaling [30]. The RLK (Bradi2gl2120.1), a 



conserved phosphoprotein among B. distachyon, O. sativa 
and A. thaliana, was phosphorylated at Serl49 (Table 2). 
Hematy et al. [65] found a RLK (THE1) mediating the 
response of growing plant cells to inhibition of cellulose 
synthesis in Arabidopsis and that it can be phosphorylated 
at a Thr residue. The phosphorylation of RLK identified 
in our study may be involved in cell wall remodeling 
during leaf development. A wall-associated receptor 
kinase (WAK, Bradi3g49160.1), as a kind of RLK, has 
the potential to serve as both linkers of the cell wall to the 
plasma membrane and as a signaling molecule. WAKs also 
function in pathogen response and cell expansion and 
elongation [27,66,67]. WAK proteins contain a cytoplasmic 
Ser/Thr kinase (STK) domain, a transmembrane domain, 
and an extra cytoplasmic region with several epidermal 
growth factor (EGF) repeats [28,68]. Interestingly, the 




Figure 7 Phosphorylation regulatory network centered by MAPK cascade. Protein kinases and protein phosphatases are highlighted in red 
and green, respectively and other phosphoproteins identified in this study are shown with sky-blue nodes. Edges with blue color represent the 
phosphoproteins with direct relationship to MEKK, MAP2K or MAPK. 



Lv et al. BMC Genomics 2014, 15:375 
http://www.biomedcentral.com/1471-2164/15/375 



Page 13 of 18 



phosphorylated site of this BdWAK was located on the C- 
terminal (Ser693) in our study, which was outside these 
domains above. As a kind of RLK, phosphorylation of this 
BdWAK may regulate functions of itself and further medi- 
ate the phosphorylation of enzymes involved in synthesis of 
cellulose in cell walls. 

The primary cell walls of grasses consist mainly of cel- 
lulose, arabinoxylan and uniquely (1, 3; 1, 4)-p-D-glucan, 
small quantities of structural proteins, and low levels of 
pectin and xyloglucan [63,69-71]. Cellulose is by far the 
most abundant biological polymer on earth, and is synthe- 
sized by cellulose synthase (CESA) contained in plasma 
membrane-localized complexes. In A. thaliana, CESA1, 
CESA3 and CESA6 interact with each other and are 
largely responsible for cellulose production during forma- 
tion of the primary cell wall [72-75], whereas three types 
of CESA subunits (CESA4, CESA7 and CESA8) are re- 
quired for secondary cell wall formation [74,76,77]. In 
addition, CESA2, CESA5, and CESA9 also function in 
primary cell wall synthesis [72]. This study identified 6 
phosphoproteins as different CESAs in B. distachyon, 
namely BdCESAl, BdCESA2, BdCESA4, BdCESA7, Bd 
CESA8 and BdCESA9 (Bradi2g34240.1, Bradilg04597.1, 
Bradi2g49912.1, Bradi3g28350.1, Bradilg54250.1 and Bradi 
4g30540.1 in Additional file 1: Table SI). Representative 
MS/MS spectra of phosphopeptides BdCESA2 and BdCE 
SA4 are shown in Additional file 10: Figure S6, panel A and 
B, respectively. AtCESA2 expressed in all organs except 
root hairs, and lack of AtCESA2 expression can lead to 
stunted growth of hypocotyls in seedlings and gready re- 
duced seed production in mature plants [78]. All cellulose 
synthases described to date have shown conserved 
structural features [78]. The amino N-terminal of the 
CESA protein contains a RING finger domain with a 
strongly conserved CxxC motif (C = Cys, x = any amino 
acid) [79]. The RING finger domain is a small zinc- 
binding domain found in many functionally distinct proteins 
[80]. Kurek et al. [81] confirmed that GhCESAl formed 
homodimers as well as heterodimers with GhCESA2 
through a zinc-binding domain. Chu et al. [78] demon- 
strated that AtCESA2 could homodimerize and the zinc 
finger domain of AtCESA2 was important for the inter- 
action. Sequence analyses indicated that all six identified 
BdCESAs contained a RING finger domain with a con- 
served CxxC motif (Additional file 10: Figure S6, panel C). 
There are 8 transmemberane regions in BdCESAl, 
BdCESA2, BdCESA8 and BdCESA9, but there is no the 
second transmemberane region in BdCESA4 and BdCE 
SA7. All phosphorylation sites identified from these 
BdCESAs are located between the RING finger domain 
and transmemberane region I. Nuhse et al. [8] performed 
a plasma membrane phosphoproteomics in Arabidopsis 
and revealed that CESA1, CESA3, and CESA5 proteins 
were phosphorylated at several sites clustered in the two 



hypervariable regions. Chen et al. [82] found that phos- 
phorylation of CESA1 differentially affects polar inter- 
action with microtubules and may regulate the length or 
quantity of a subset of cellulose microfibrils in the primary 
cell wall. The region with phosphorylation sites in each 
BdCESA belonged to hypervariable region I (Additional 
file 10: Figure S6, panel C). Phosphorylation of AtCESA7 
was linked to its degradation via a 26S proteasome- 
dependent pathway [83]. It seems that phosphorylation of 
the 6 BdCESAs may play vital roles in cell expansion dur- 
ing rapid vegetative growth of Brachypodium. 

(1, 3; 1, 4)-|3-D-glucan is mainly detected in some fam- 
ilies of the Poales [70]. Recently, various genes involved in 
(1, 3; 1, 4)-[3-D-glucan synthesis were identified [84,85]. 
They belong to two grass-specific families (CSLF and 
CSLH) of the cellulose synthase-like (CSL) gene super- 
family. In the present study, two CSL proteins BdCSLF6 
(Bradi3gl6307.1) and BdCSLHl (Bradi5gl0130.1) were 
identified as phosphoproteins. In developing barley 
coleoptiles, the transcript levels of HvCslF6 are maximal 
in 4- to 5-d coleoptiles, at the time when (1, 3; 1, 4)-|3- 
D-glucan content of coleoptile cell walls also reaches its 
maximum level [84]. HvCSLHl protein is responsible 
for p-glucan deposited in cell walls of transgenic Arabi- 
dopsis [85]. In B. distachyon, the CSLH gene showed a 
much higher transcript level and played a prominent 
role in (1, 3; 1, 4)-p-D-glucan synthesis [64]. 

Conclusions 

In this work, TiGVenriched phosphopeptides in Bd21 
seedling leaves detected by LC-MS/MS, and MaxQuant 
software were used to identify the phosphoproteins and 
phosphorylation sites. To our knowledge, this is the first 
large-scale phosphoproteome analysis in B. distachyon and 
provides an overview of in vivo phosphorylation events in 
seedling leaves. A total of 950 phosphoproteins containing 
1470 unambiguous phosphorylation sites were identified. 
In-depth GO, conservation, and KEGG pathway analyses 
of these phosphoproteins revealed the phosphorylation 
profiling during rapid seedling leaf growth. GO analysis 
showed that the 950 phosphoproteins mainly functioned 
in protein binding, nucleotide binding and kinase activity. 
Conservation analysis revealed that 887 (93.4%) of the 
phosphoproteins had phosphorylation orthologs in other 
species. Many fundamental biological pathways and two 
signal transduction pathways were identified through 
KEGG pathway analysis. 

Fifty eight TFs, 94 PKs and 8 PPs were identified from 
the large dataset and analyzed in depth. Thirteen enriched 
motifs specific for various Ser/Thr PKs were identified 
from B. distachyon. A total of 127 phosphoproteins con- 
tained multiple phosphorylation sites and 91 of them 
contained two or more phosphorylation motifs. Several 
signal transduction and cell wall expansion associated 
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phosphoproteins that may have vital functions in the 
rapid growth of B. distachyon seedling leaves were 
found such as B-ARR, PTEN, CESA and CSL. Our over- 
all results provide an important resource for further 
research and new insights into the regulation of protein 
post-translational phosphorylation during rapid vegeta- 
tive growth of Poaceae species. 

Methods 

Plant materials 

Seeds of Brachypodium distachyon Bd21, kindly pro- 
vided by Dr. John Vogel, USDA-ARS, Albany, CA, were 
surface sterilized in 5% sodium hypochlorite for 5 min, 
and rinsed 4 times in sterile distilled water. Seeds were 
submerged in water for 12 h at room temperature, and 
then transferred to wet filter paper for 24 h to germinate 
at room temperature (22-25°C). Uniformly germinated 
seeds were selected and grown in three plastic pots con- 
taining Hoagland's solution [86]; the Hoagland's solution 
was changed every two days. At the three leaf stage, all 
three leaves were collected and frozen at -80°C. Culture 
of the Bd21 seedlings was repeated three times and the 
leaves were collected independently. 

Protein extraction 

Total proteins were extracted from seedling leafs accord- 
ing to the method of Wang et al. [87] with minor modi- 
fications. Approximately 400 mg fresh leaves of each 
sample were ground into fine powder in liquid nitrogen. 
The ground powder was suspended in 4 mL SDS buffer 
(30% sucrose, 2% SDS, 100 mM Tris-HCl, pH 8.0, 50 mM 
EDTA-Na 2 , 20 mM DTT) and 4 mL phenol (Tris-buff- 
ered, pH 8.0) in a 10 mL tube, and 1 mM phenylmetha- 
nesulfonyl fluoride (PMSF) and PhosSTOP Phosphatase 
Inhibitor Cocktail (Roche, Basel, Switzerland), were added 
to inhibit the activities of proteases and phosphatases. The 
mixtures were thoroughly vortexed thoroughly for 30 s 
and the phenol phase was separated by centrifugation at 
14,000 xg and at 4°C for 15 min. The upper phenol phase 
was pipetted to fresh 10 mL tubes and four volumes of 
cold methanol plus 100 mM ammonium acetate were 
added, and the mixture was stored at -20°C for at least 
30 min. After centrifugation at 14,000 x g and at 4°C for 
15 min, the supernatant was carefully discarded and the 
precipitated proteins were washed twice with cold metha- 
nolic ammonium acetate (100 mM) and ice-cold 80% acet- 
one, respectively. Finally the pellet was vacuum-dried and 
then dissolved in lysis buffer (7 M urea, 2 M thiourea, 
4% w/v CHAPS and 65 mM DTT) over 3 h at 4°C The 
protein mixtures were harvested by centrifugation at 
14,000 x g and 4°C for 15 min to remove insoluble mate- 
rials. The concentrations of the extracted protein mix- 
tures were determined with a 2-D Quant Kit (Amersham 
Bioscience, Buckinghamshire, UK) using BSA (2 mg/mL) 



as standard, and the final protein solution was stored at - 
80°C for later use. 

Phosphopeptide enrichment using Ti0 2 microcolumns 

Extracted protein mixtures were directly reduced with 
dithiothreitol (DTT), alkylated with iodoacetamide, and 
subsequently digested with endoproteinase Lys-C and 
trypsin as previously described [88]. The enrichment pro- 
cedure for phosphopeptides was performed as reported 
by Wu et al. [5] with modifications. The TiC>2 beads (GL 
Sciences, Tokyo, Japan) were incubated in 400 fiL loading 
buffer containing 65% Acetonitrile (ACN)/2% trifluoroace- 
ticacid (TFA)/saturated by glutamic acid. A total of 2 mg 
of tryptic peptides were dissolved in 600 [A. loading buffer, 
and incubated with an appropriate amount of TiC>2 beads. 
After washing with 600 fih wash buffer (65% ACN/0.1% 
TFA), phosphopeptides were eluted twice with 300 ^L elu- 
tion buffer (500 mM NH4OH/60% ACN) and the eluates 
were dried down and reconstituted in 0.1% formic acid 
(FA)/H 2 0 for MS analysis. 

Phosphopeptide identification and phosphorylation site 
localization using LC-MS/MS 

Enriched phosphopeptides were separated on a self- 
packed C18 reverse phase column (75 fim I.D., 150 mm 
length) (Column Technology Inc., Fremont, CA), which 
directly connected the nano electrospray ion source to a 
LTQ-Orbitrap XL mass spectrometer (Thermo Fisher 
Scientific, San Jose, CA). Pump flow was split to achieve 
a flow rate at 1 ^L/min for sample loading and 300 nL/ 
min for MS analysis. The mobile phases consisted of 
0.1% FA (A), and 0.1% FA and 80% ACN (B). A five-step 
linear gradient of 5% to 30% B in 105 min, 35% to 90% B 
in 16 min, 90% B in 4 min, 90% to 2% B for 0.5 min and 
2% B for 14.5 min was performed. The spray voltage was 
set to 2.0 kV and the temperature of the heated capillary 
was 240°C. 

For data acquisition, each MS scan was acquired at a 
resolution of 60,000 (at 400 m/z) with lock mass option 
enabled and was followed by a data-dependent top 10 
MS/MS scans using collision induced dissociation (CID). 
The threshold for precursor ion selection was 500 and 
mass window for precursor ion selection was 2.0 Da. 
The dynamic exclusion duration was 120 s, repeat count 
was 1 and repeat duration was 30 s. The analyzer for the 
MS scans was Orbitrap and for the MS/MS scans LTQ 
(37% relative collision energy). Three biological repli- 
cates were performed independently from sample collec- 
tion to phosphopeptide identification using LC-MS/MS. 

The raw files were processed with MaxQuant (version 
1.2.2.5) [89] and searched against the B. distachyon protein 
database (31,029 entries in total) in Phytozome (http:// 
www.phytozome.net/search.php; version 9.1) concatenated 
with a decoy of reversed sequences. The following 
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parameters were used for database searches: cysteine 
carbamidomethylation was selected as a fixed modifica- 
tion; methionine oxidation, protein N-terminal acetylation, 
and phosphorylation on serine, threonine and tyrosine 
were selected as variable modifications. Up to two missing 
cleavage points were allowed. The precursor ion mass tol- 
erances were 7 ppm, and fragment ion mass tolerance was 
0.5 Da for MS/MS spectra. The false discovery rate (FDR) 
was set to < 1.0% for both peptide and protein identifica- 
tions, the minimum peptide length was set to 6. 

Phosphorylation site localization was based on PTM 
scores that assign probabilities for each of the possible 
sites according to their site-determining ions. In this 
study, MaxQuant (version 1.2.2.5) was used to calculate 
PTM scores and PTM localization probabilities. Poten- 
tial phosphorylation sites were then grouped into three 
categories depending on their PTM localization probabil- 
ities [10,37,88], namely class I (localization probability, 
P > 0.75), class II (0.75 > P > 0.5) and class III (P < 0.5). A 
false discovery rate (FDR) of 1% was used for phosphoryl- 
ation sites identification. 

Bioinformatics 

The GO biological processes, molecular functions and 
cellular components of the identified phosphoproteins 
were examined by using Blast2GO software [90]. GO 
enrichment analysis was conducted with agriGO [91]. 
The statistical test method was set as Fisher and the 
multi-test adjustment method was set as Bonferroni. 
Pfam domain information was extracted from the data- 
base (http://pfam.sanger.ac.uk/search) [92]. Protein de- 
scriptions were extracted from the B. distachyon protein 
database (31,029 entries in total) in Phytozome (http:// 
www.phytozome.net/search.php; version 9.1). For path- 
way analysis, the proteins were searched in the KEGG 
B. distachyon database and mapped to B. distachyon- 
specific pathways with KEGG Mapper. The enrichment of 
TF family significantly regulated by phosphorylation was 
performed using a two-sided test of the hypergeometric 
distribution. The significantly enriched phosphorylation 
motifs were extracted from phosphopeptides with unam- 
biguous phosphorylation sites (class I) using the Motif-X 
algorithm [33]. The phosphopeptides were centered at 
phosphorylated amino acid residues and aligned, and six 
positions upstream and downstream of the phosphoryl- 
ation site were included. In the cases of C- and N- 
terminal peptides, the sequence was completed to 13 
amino acids with the required number of "X"s, where X 
represents any amino acid. Because the upload restric- 
tion of Motif-X is 10 MB, a fasta format dataset (nearly 
10 MB) containing the protein sequences from the B. 
distachyon protein database in Phytozome (version 9.1) 
was used as the background database to normalize 
scores against random distributions of amino acids. The 



occurrence threshold was set to 5% of the input data set at 
a minimum of 20 peptides, and the probability threshold 
was set to P < 10~ . The sequences of all the phospho- 
proteins were used for BLAST analysis with the Na- 
tional Center for Biotechnology Information (NCBI) 
clusters of the KOG database to obtain the KOG numbers 
of those proteins. Then a dataset containing all the KOG 
numbers was used for protein-protein interaction (PPI) 
analysis by using the STRING database (version 9.1, 
http://string-db.org) [93]. Only interactions that had a 
confidence score of at least 0.9 were used to construct the 
network and the networks were then displayed using 
Cytoscape (version 3.0.2) software [94]. HMMTOP (ver- 
sion 2.0, http://www.enzim.hu/hmmtop/) was used to pre- 
dict the transmembrane helices in proteins [95]. 

Additional files 



Additional file 1: Table SI. Phosphorylation sites and phosphoproteins 
identified from seedling leaves of B. distachyon. (A) Total phosphorylation 
sites and phosphoproteins. (B) Annotation of the phosphoproteins. 

Additional file 2: Figure SI. Representative tandem mass spectrometry 
spectra of the seven phosphopeptides containing the eight phosphosites 
in Bradilg66870.1. (A) KHLPSQDITDAS(ph)GDEEEGSR. (B) HWS(ph)PSAGRR. 
(C) LKDFS(ph)ADPELNK. (D) MNYLGT(ph)PPS(ph)DLEK. (E) NS(ph)PISIK. (F) 
QFDS(ph)PDDSLVK. (G) RIS(ph)PPSGR. 

Additional file 3: Table S2. Conservation analysis of phosphoproteins 
identified from B. distachyon. 

Additional file 4: Figure S2. GO biological process enrichment of 582 
highly conserved phosphoproteins. The statistical significance of the 
enrichment analysis is represented by a scale of red tones whose intensity is 
proportional to the degree of significance starting from FDR < 0.05. 

Additional file 5: Table S3. Pathway analysis of the identified 
phosphoproteins by KEGG. (A) Pathway search results of all identified 
phosphoproteins and highly conserved phosphoproteins. (B) KEGG 
identifiers of all the identified phosphoproteins. 

Additional file 6: Figure S3. Pathway mapping of identified 
phosphoproteins by KEGG. Phosphoproteins were analyzed in the KEGG 
B. distachyon pathway database. The phosphoproteins were used to 
search the KEGG B. distachyon database and mapped to B. distachyon 
specific pathways with KEGG Mapper. In each pathway map, objects with 
red foreground color and blue background color represent the highly 
conserved phosphoproteins; objects with pink foreground and blue 
background represent the other conserved phosphoproteins; objects 
with yellow foreground and blue background represent proteins 
identified as novel phosphoproteins; objects with green background 
represent the total proteins in the S. distachyon KEGG database. Details of 
the pathways can be found in Additional file 5: Table S3. 

Additional file 7: Table S4. Enriched phosphorylation motifs, and the 
phosphopeptides and the corresponding phosphoproteins in each 
phosphorylation motif. (A) Phosphorylation motifs enriched by Motif-X and 
the putative protein kinases. (B) Phosphopeptides and the corresponding 
phosphoproteins in each phosphorylation motif. (C) Motifs of 
phosphoproteins with multiple phosphorylation sites. 

Additional file 8: Figure S4. GO enrichment of phosphoproteins 
represented by the phosphorylation motifs in B. distachyon. GO 
enrichment graphs of phosphoproteins contained each of 1 1 motifs were 
displayed by Blast2GO software. No GO terms were significantly enriched 
from phosphoproteins containing motif 7 or motif 1 1. The statistical 
significance of the enrichment analysis is represented by a scale of red 
tones whose intensity is proportional to the degree of significance 
starting from p < 0.05. 
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Additional file 9: Figure S5. Protein-protein interaction network of the 
phosphoproteins. The confidence (score) was set as the highest (0.900) 
and only the phosphoproteins identified in this study were used to 
construct the interaction profiles. Protein kinases and protein phosphatases 
are highlighted in red and green respectively and other phosphoproteins 
identified in this study are shown with sky-blue nodes. 

Additional file 10: Figure S6. Phosphorylation sites of the six BdCESAs. 
(A) and (B) The tandem mass spectrometry spectra of the phosphopeptides 
EFSGS(ph)LGNVAWK from BdCESA2 (Bradilg04597.1) and VTIASQLS(ph)DR 
from BdCESA4 (Bradi2g4991 2.1 ). (C) Seguence alignment among the six 
BdCESAs. Phosphorylation sites identified in this study are highlighted in 
red. The RING domain and transmembrane regions are highlighted in 
yellow and green, respectively. 
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